cd "C:\insert file"

local now = subinstr("`c(current_time)'", ":", "", .)
display `now'


/*
This code compares three strategies for prospectively assigning participants to conditions:
1. random assignment
2. stratification (i.e., using separate lists of counter-balanced random numbers to assign participants to conditions)
3. minimization procedures
*/

/*
Users should:
1. adjust the "user defined variables below, as needed
2. run the code, which will produce one dataset for each type of assignment, named:
	a. results_randomAssignment
	b. results_stratify
	c. results_minimization
*/

//*******************************************
//*******************************************
//USER-DEFINED VARIABLES
local numOfExperimentalConditions = 16
local sampleSizeOfInterest = 100
local numberOfSimulations = 100
//NOTE: stratify variables must be binary
local stratifyVar1 = "d_langSpanish"
local stratifyVar2 = "d_langHaitian"
local stratifyVar3 = "d_medicaid"
local stratifyVar4 = "d_ASD"
local minimPercentRandom = .1 //% of assignments made at random
local minimVar1 = "total" //always keep this variable--balances total size of each condition
local minimVar2 = "d_medicaid"
local minimVar3 = "d_workOutsideHome"
local minimVar4 = "d_female"
local minimVar5 = "d_langSpanish"
local minimVar6 = "d_langHaitian"
local minimVar7 = "ageref"
local minimVar8 = "d_ASD"
local minimVar9 = "d_HSgrad"
local minimVar10 = "d_molp"


//*******************************************
//location of saved results

mkdir "Sim`numberOfSimulations' n`sampleSizeOfInterest' `c(current_date)' `now'"
cd "Sim`numberOfSimulations' n`sampleSizeOfInterest' `c(current_date)' `now'"

//log using "C:\file subinstr("`c(current_time)'", ":", "", .) `c(current_date)'.smcl"
log using "Log `c(current_date)' `now'.smcl",replace
display "simulation began on `c(current_date)' at `c(current_time)'"

//*******************************************
//*******************************************

	//declare matrices to save summary results
	matrix rMinRnd = J(`numberOfSimulations',10,.)
	matrix rMaxRnd = J(`numberOfSimulations',10,.)
	matrix rRangeRnd = J(`numberOfSimulations',10,.)
	matrix rMinStrat = J(`numberOfSimulations',10,.)
	matrix rMaxStrat = J(`numberOfSimulations',10,.)
	matrix rRangeStrat = J(`numberOfSimulations',10,.)
	matrix rMinMnmztn = J(`numberOfSimulations',10,.)
	matrix rMaxMnmztn = J(`numberOfSimulations',10,.)
	matrix rRangeMnmztn = J(`numberOfSimulations',10,.)

	//open master data
	clear
	import delimited "C:File.txt"
	gen d_medicaid = 0
	replace d_medicaid = 1 if insurance=="Public Insurance (Medicaid)"
	
	gen d_workOutsideHome = 0 
	replace d_workOutsideHome = 1 if work_outside_home=="Yes"
	
	gen d_female = 0 
	replace d_female = 1 if gender=="Female"

	gen d_langSpanish=0
	replace d_langSpanish=1 if referral_language=="Spanish"

	gen d_langHaitian=0
	replace d_langHaitian=1 if referral_language=="Haitian Creole"
	
	gen d_ASD=0
	replace d_ASD=1 if diagnosis=="ASD"

	gen d_HSgrad = 0 
	replace d_HSgrad = 1 if hsgrad=="Yes"

	gen d_molp = 0 
	replace d_molp = 1 if molp=="Yes"

	gen id = reference_id
	
	gen total = 1
	//store maximum value for each minimization variable:
	//	1 = binary
	//	2 = ordinal
	//	0 = no variable declared
	local numberOfMinimizationVars = 0 
	forvalues j = 1/10 {
		quietly sum `minimVar`j''
		local mean = r(mean)
		local minimVar`j'Max = cond("`minimVar`j''" != "", r(max),0)
		display `minimVar`j'Max'
		//if not binary, create binary variable coded at the mean
		//NOTE: there are other ways to address non-binary variables in minimization that we can explore
		if `minimVar`j'Max' !=0 & `minimVar`j'Max' !=1 {
			gen `minimVar`j''B = `minimVar`j''
			local rx = `mean'+.00000001
			recode `minimVar`j'' (min/`mean' = 0) (`rx'/max = 1)
			//local minimVar`j' = "minimVar`j'B"
			}
		local numberOfMinimizationVars = cond("`minimVar`j''" != "",`numberOfMinimizationVars'+1,`numberOfMinimizationVars')
		}
		


	save "C:\file.dta", replace


forvalues i = 1/`numberOfSimulations' {

	//PRELIMINARY CODE
	clear
	use "C:\file.dta"
	gen order = .
	/*
	//create variable of random numbers and order by that variable
	//"order" represents the order in which participants present to the study
	gen order = runiform()
	sort order

	//keep number of participants corresponding to "sampleSizeOfInterest"
	drop if _n > `sampleSizeOfInterest'

	//save "C:\file - run `numberOfSimulations'", replace
	*/
	//take bootstrap sample of desired sample size
	bsample `sampleSizeOfInterest'

	//*******************************************
	//SIMULATE RANDOM ASSIGNMENT
	//*******************************************
	gen expCondition_randomAssignment = trunc(runiform(0,`numOfExperimentalConditions')+1)
	/*
	tab expCondition_randomAssignment
	*/
	
	//*******************************************
	//SIMULATE STRATIFICATION
	//*******************************************	
	//	This code defines the number of strata by the stratification variables assigned above (2^#vars)
	//	Within each stratum, the participants are sorted by a random number
	//	assignments to conditions are then made sequentially, beginning with a random number between 1 and the number of conditions
	
	//	This process is designed to simulated use of lists of random numbers that are chosen to ensure equal sample sizes at a given "sampleSizeOfInterest"
	//  Note that in reality, use of balanced lists of random numbers may slightly underperform the method simulated here
	
	egen strata = group(`stratifyVar1' `stratifyVar2' `stratifyVar3' `stratifyVar4')
	sum strata
	local numStrata = r(max)

	gen expCondition_stratify = .
	
	forvalues j = 1/`numStrata' {
		local begin = runiformint(1,`numOfExperimentalConditions')
		gen temp_random = runiform(0,1) if strata == `j'
		sort temp_random
		sum id if strata == `j'
		local sizeOfStratum = r(N)
		
		local s = 1
		local n = `begin'
		forvalues k = 1/`sizeOfStratum' {
		display "n: " `n'
		display "s: " `s'
			replace expCondition_stratify = `n' if _n == `k'
			local s = `s' + 1
			local n = `n' + 1
			local s = cond(`s'>`numOfExperimentalConditions',1,`s')
			local n = cond(`n'>`numOfExperimentalConditions',1,`n')
			local n = cond(`s'>`numOfExperimentalConditions',`begin',`n')
			}
		drop temp_random
		}
	sort order

	//*******************************************
	//SIMULATE MINIMIZATION
	//*******************************************	

	//cycle through each participant one at a time
		matrix tempResults = J(`numOfExperimentalConditions'+1,11,0)
		matrix colnames tempResults = `minimVar1' `minimVar2' `minimVar3' `minimVar4' `minimVar5' `minimVar6' `minimVar7' `minimVar8' `minimVar9' `minimVar10' rangeSum
		

	gen expCondition_minimize = .
	forvalues j = 1/`sampleSizeOfInterest' {
		display "simulation `i' of `numberOfSimulations'"
		display "participant `j' of `sampleSizeOfInterest'"
		//1. consider assignment to each condition, 
		//2. calculate range (i.e., max frequency - min frequency) across conditions with respect to each minimization variable;
		//3. sum values of range across variables
		//4. indicate assignment to condition that yields smallest difference
		//5. in the event of a tie, select at random
		local expConditionIndicated = 1
		local minSumDiff = 1000000000
		quietly {
				forvalues k = 1/`numOfExperimentalConditions' {
					//replace expCondition_minimize = `k' if _n == `j'
					local sum = 0 
					local count = 0 
					//for each minimization variable...
					forvalues m = 1/`numberOfMinimizationVars' {
						//calculate frequency of positive value for each variable in each condition (used later to calculate range)
						forvalues n = 1/`numOfExperimentalConditions' {
							capture sum `minimVar`m'' if expCondition_minimize==`n'
							local freq = r(mean)*r(N)
							local freq = cond(`freq'==.,0,`freq')
							matrix tempResults[`n',`m'] = `freq'
							//local sum = `sum' + `mean'
							//local count = cond("`minimVar`m''" != "",`count'+1,`count')
							}
						//matrix tempResults[`numOfExperimentalConditions'+1,`m'] = `sum'/`count'
						}
					
					//1. consider assignment to each condition...	
					local sum = 0 
					forvalues q = 1/`numOfExperimentalConditions' {
						matrix tmp = tempResults
						/*
						noisily display "minimvar1: " `minimVar1'[`j']
						noisily display "minimvar2: " `minimVar2'[`j']
						noisily display "minimvar3: " `minimVar3'[`j']
						noisily display "minimvar4: " `minimVar4'[`j']
						noisily display "minimvar5: " `minimVar5'[`j']
						noisily display "minimvar6: " `minimVar6'[`j']
						noisily display "minimvar7: " `minimVar7'[`j']
						noisily display "minimvar8: " `minimVar8'[`j']
						*/

						matrix tmp[`q',1] = cond(`minimVar1'[`j']==1,tmp[`q',1] + 1,tmp[`q',1])
						matrix tmp[`q',2] = cond(`minimVar2'[`j']==1,tmp[`q',2] + 1,tmp[`q',2])
						matrix tmp[`q',3] = cond(`minimVar3'[`j']==1,tmp[`q',3] + 1,tmp[`q',3])
						matrix tmp[`q',4] = cond(`minimVar4'[`j']==1,tmp[`q',4] + 1,tmp[`q',4])
						matrix tmp[`q',5] = cond(`minimVar5'[`j']==1,tmp[`q',5] + 1,tmp[`q',5])
						matrix tmp[`q',6] = cond(`minimVar6'[`j']==1,tmp[`q',6] + 1,tmp[`q',6])
						matrix tmp[`q',7] = cond(`minimVar7'[`j']==1,tmp[`q',7] + 1,tmp[`q',7])
						matrix tmp[`q',8] = cond(`minimVar8'[`j']==1,tmp[`q',8] + 1,tmp[`q',8])
						
						//noisily display "tmp[`q',9]: " tmp[`q',9]
						//noisily display "minimvar9: " `minimVar9'[`j']
						matrix tmp[`q',9] = cond(`minimVar9'[`j']==1,tmp[`q',9] + 1,tmp[`q',9])
						//noisily display "tmp[`q',9]: " tmp[`q',9]

						matrix tmp[`q',10] = cond(`minimVar10'[`j']==1,tmp[`q',10] + 1,tmp[`q',10])
						local sum = 0
						forvalues m = 1/`numberOfMinimizationVars' {
							local min = `sampleSizeOfInterest' 
							local max = 0 
							
							forvalues n = 1/`numOfExperimentalConditions' {
								local max = cond(tmp[`n',`m'] > `max',tmp[`n',`m'],`max')
								local min = cond(tmp[`n',`m'] < `min',tmp[`n',`m'],`min')
								}
							
							local range = `max'-`min'
							matrix tmp[`numOfExperimentalConditions'+1,`m'] = `range'
							local sum = `sum' + (`range')
							}
							matrix tempResults[`q',11] = `sum'
						}
					}
						//identify and count conditions with minimum sum
						//identify minimum
						local min = .
						forvalues n = 1/`numOfExperimentalConditions' {
								local min = cond(tempResults[`n',11] < `min',tempResults[`n',11],`min')
								}	
						
							noisily matrix list tempResults
							noisily display "The minimum sum of ranges is: " `min' " [NOTE: the final row is not a condition]"
						//count conditions with minimum
						local count = 0 
						forvalues n = 1/`numOfExperimentalConditions' {
								local count = cond(tempResults[`n',11] == `min',`count'+1,`count')
								}	
							noisily display "There are " `count' " conditions with this minimum value"						
						//randomly select among these conditions and assign case to condition
						local rand = runiformint(1, `count')
						local count2 = 0 
						forvalues n = 1/`numOfExperimentalConditions' {
								local count2 = cond(tempResults[`n',11] == `min',`count2' + 1,`count2')
									//noisily display "count2: " `count2'
								local expConditionIndicated = cond(`count2'==`rand' & tempResults[`n',11] == `min',`n',`expConditionIndicated')	
								}	
					
					
				}
			noisily display "The participant is assigned at random to the " `rand' "th one..."
			noisily display "Which is condition number " `expConditionIndicated'
		replace expCondition_minimize = `expConditionIndicated'  if _n == `j'
		//matrix list tempResults

		}

	//save final results to matrices for random assignment	
	tabstat `minimVar1' `minimVar2' `minimVar3' `minimVar4' `minimVar5' `minimVar6' `minimVar7' `minimVar8' `minimVar9' `minimVar10', by(expCondition_randomAssignment) stats(sum) save
	forvalues k = 1/`numberOfMinimizationVars' {
		local min = `sampleSizeOfInterest'
		local max = 0
		forvalues m = 1/`numOfExperimentalConditions' {	
			matrix tmp = r(Stat`m')
			local min = cond(tmp[1,`k'] < `min', tmp[1,`k'], `min')
			local max = cond(tmp[1,`k'] > `max' & tmp[1,`k'] !=., tmp[1,`k'], `max')
			}
		local range = `max' - `min'
		matrix rMinRnd[`i',`k'] = `min'
		matrix rMaxRnd[`i',`k'] = `max'
		matrix rRangeRnd[`i',`k'] = `range'
		}
		
	//save final results to matrices for stratification	
	tabstat `minimVar1' `minimVar2' `minimVar3' `minimVar4' `minimVar5' `minimVar6' `minimVar7' `minimVar8' `minimVar9' `minimVar10', by(expCondition_stratify) stats(sum) save
	forvalues k = 1/`numberOfMinimizationVars' {
		local min = `sampleSizeOfInterest'
		local max = 0
		forvalues m = 1/`numOfExperimentalConditions' {	
			matrix tmp = r(Stat`m')
			local min = cond(tmp[1,`k'] < `min', tmp[1,`k'], `min')
			local max = cond(tmp[1,`k'] > `max' & tmp[1,`k'] !=., tmp[1,`k'], `max')
			}
		local range = `max' - `min'
		matrix rMinStrat[`i',`k'] = `min'
		matrix rMaxStrat[`i',`k'] = `max'
		matrix rRangeStrat[`i',`k'] = `range'
		}
		
	//save final results to matrices for minimization	
	tabstat `minimVar1' `minimVar2' `minimVar3' `minimVar4' `minimVar5' `minimVar6' `minimVar7' `minimVar8' `minimVar9' `minimVar10', by(expCondition_minimize) stats(sum) save
	forvalues k = 1/`numberOfMinimizationVars' {
		local min = `sampleSizeOfInterest'
		local max = 0
		forvalues m = 1/`numOfExperimentalConditions' {	
			matrix tmp = r(Stat`m')
			local min = cond(tmp[1,`k'] < `min', tmp[1,`k'], `min')
			local max = cond(tmp[1,`k'] > `max' & tmp[1,`k'] !=., tmp[1,`k'], `max')
			}
		local range = `max' - `min'
		matrix rMinMnmztn[`i',`k'] = `min'
		matrix rMaxMnmztn[`i',`k'] = `max'
		matrix rRangeMnmztn[`i',`k'] = `range'
		}
	}


local filename = "resultsFromLog `c(current_date)' `now'.xlsx"
putexcel set "`filename'", sheet(rMinRnd) modify
putexcel A1 = matrix(rMinRnd)
putexcel set "`filename'", sheet(rMaxRnd) modify
putexcel A1 = matrix(rMaxRnd)
putexcel set "`filename'", sheet(rRangeRnd) modify
putexcel A1 = matrix(rRangeRnd)
putexcel set "`filename'", sheet(rMinStrat) modify
putexcel A1 = matrix(rMinStrat)
putexcel set "`filename'", sheet(rMaxStrat) modify
putexcel A1 = matrix(rMaxStrat)
putexcel set "`filename'", sheet(rRangeStrat) modify
putexcel A1 = matrix(rRangeStrat)
putexcel set "`filename'", sheet(rMinMnmztn) modify
putexcel A1 = matrix(rMinMnmztn)
putexcel set "`filename'", sheet(rMaxMnmztn) modify
putexcel A1 = matrix(rMaxMnmztn)
putexcel set "`filename'", sheet(rRangeMnmztn) modify
putexcel A1 = matrix(rRangeMnmztn)

matrix list rRangeRnd	
matrix list rRangeStrat	
matrix list rRangeMnmztn	

display "simulation completed on `c(current_date)' at `c(current_time)'"

log close